# TUTORIAL 5 : CO2 adsorption in a zeolite¶

Authors: James Grant (r.j.grant{at}bath.ac.uk), Tom L. Underwood (t.l.underwood{at}bath.ac.uk)

## Introduction¶

In this exercise you will evaluate the amount of CO2 adsorbed in a siliceous zeolite under certain conditions. We use a Grand Canonical Monte Carlo simulation, in which Monte Carlo moves are performed which insert and remove CO2 molecules from the system, alongside translational moves to explore the local environment. After a certain amount of simulation time, the number of CO2 molecules will reflect the equilibrium concentration within the zeolite. You will examine the length of time required to run simulations under different conditions and how to analyse the data in order to construct an adsorption isotherm. Of course, there will be fluctuations in the number of CO2 molecules about the ‘true’ average which you can also investigate.

## CONFIG¶

As always the CONFIG file contains the starting structure. If you examine the file you will see that the system contains two ‘molecule’ types: ‘zeolite’, which contains silcon and oxygen atoms (labeled Si and O) making up the zeolite and Xe atoms (labelled Xe). and ‘co2’, which contains C (labeled C_) and O (labeled O_). There is one ‘zeolite’ molecule and 52 CO2 molecules. Your first objective is to visualise the structure and identify why the zeolite contains Xe atoms. HINT: Also check the FIELD file…

## FIELD¶

The FIELD file is shown below:

  1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 Force fields and bond constraints for for CO2 in a zeolite (use EPM2 for VDW) CUTOFF 12.0 UNITS kcal NCONFIGS 1 ATOMS 5 Si core 28 2.4 O_ core 16 -1.2 C_ core 12 0.6512 O_C core 16 -0.3256 Xe core 1 0.0 MOLTYPES 2 zeolite MAXATOM 584 co2 ATOMS 3 3 C_ core 0.00000000 0.0000000 0.0000000 O_C core 1.16300000 0.0000000 0.0000000 O_C core -1.16300000 0.0000000 0.0000000 BONDS 3 1 2 1 1 3 1 2 3 2 FINISH VDW 6 O_ core C_ core lj 0.094566874 2.892 O_ core O_C core lj 0.15998351 3.033 C_ core C_ core lj 0.055892323 2.757 C_ core O_C core lj 0.094566874 2.892 O_C core O_C core lj 0.15998351 3.033 O_C core Xe core 12-6 16777216 0.0 BONDS 2 C_ core O_C core buck 0.0 0.1 0.0 O_C core O_C core buck 0.0 0.1 0.0 CLOSE 

The first thing to note is that this system has charged atoms: Si and O_ in the zeolite, C and O in CO2, and Xe, are assigned net charges of 2.4, -1.2, 0.72, -0.36 and 0 respectively. The Si and O_ potential is based on the BKS [1] potential. Formally this uses a Buckingham potential but in this instance, for the purpose of generating cross terms for interaction with the CO2 model, this has been fitted with a Lennard-Jones potential. The CO2 is modelled with the EPM [2] model. The final vdw interaction is between the CO2 oxygen atoms and the ‘xenon’ particles. These only have a strong repulsive interaction with the CO2 oxygens. The ‘xenon’ particles prevent CO2 molecules being inserted by into regions which are inaccessible by diffusion.

## CONTROL¶

The CONTROL file is shown below:

  1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 GCMC simulation of CO2 in zeolite use gaspressure # input is in partial pressure in GCMC moves # (as opposed to chemical potential) finish temperature 273.0 acceptmolmoveupdate 200 # Period (in moves) at which the # maximum move size is recalculated acceptmolrotupdate 200 # Period (in moves) at which the # maximum rotation angle is updated steps 100000 # Number of moves to perform in simulation equilibration 50000 # Equilibration time before statistics # are gathered (in moves) print 1000 # Information is output every 'print' moves revconformat dlmonte # REVCON file is in DL_POLY CONFIG format archiveformat dlpoly4 # ARCHIVE/trajectory file format # NB. No support yet for GCMC tracjectories stack 10000 # Size of blocks (in moves) for block averaging maxmolrot 0.005 # Initial maximum rotation angle (degrees) move molecule 1 20 # Perform translation moves for # 1 molecule type weight 20 co2 # move co2 molecule type move rotatemol 1 20 # Perform rotation moves for # 1 molecule type weight 20 co2 # rotate co2 molecule type move gcinsertmol 1 60 0.5 # Perform insertion/removal moves # with 1 molecule type weight 60 # with a min. distance of 0.5 from existing # atoms for inserts co2 0.0001 # Insert/remove CO2 with a # partial pressure of 0.0001 katm start 

Note that we are now using molecular moves which translate (move molecule) and rotate (move rotatemol) the CO2 molecules, but we do not use atomic moves which translate the atoms within the molecules. The CO2 are therefore considered as rigid molecules during the simulation. This restriction typically has to be in place for standard GCMC in order to satisfy detailed balance. GC moves are introduced by the move gcinsertmol directive. As with the other move types the first number is the number of molecules we will be inserting/deleting and the second is the weight given to GC moves. An additional feature is the ability to automatically reject moves with would result in any atom of the molecule lying within a specified distance of any atom already in the system. In order to satisfy detailed balance this should be smaller than the closest distance you would expect atoms to approach.

## Partial Pressure, Activity and Chemical Potential¶

The use gaspressure directive specified at the beginning of the CONTROL file means that the partial pressure of the gas, rather than the activity is specified.

$a = \gamma \frac{P}{P_0},$

where a is the activity, $$\gamma=1$$ is the gas fraction assumed to be 1, and $$P$$ and $$P_0$$ the pressure and reference pressure respectively. Note that DL_MONTE introduces a further converesion factor, dividing by 0.163882576 to convert to internal units.

The activity relates to chemical potential according to

$a = \exp((\mu - \mu_0)/RT)$

where $$\mu$$ and $$\mu_0$$ are the chemical potential and refernce chemical potential, R gas constant and T temperature.

## Exercise: Trace an Adsorption Isotherm¶

Run the intial inputs by launching DL_MONTE. We are interested in the number of adsorbed CO2 molecules. The average number and fluctuations are detailed at the end of the OUTPUT.000 and you should take note of these. However you can also extract a time sequence of molecule numbers by running the script:

strip_adsorb.sh


The results of this time sequence are output to the file ‘adsorb.dat’. At the conditions specified in the CONTROL file the number of adsorbed molecules is small, because this corresponds to a relatively low partial pressure of CO2. As we increase the partial pressure the number of CO2 molecules adsorbed will increase. Making sure you copy your output to different folders for each calculation. Increase the partial pressure until the zeolite is saturated. Compare the time sequences and final averages and fluctuations to make an initial estimate of the adsorption isotherm.

What do these initial results show would you need to consider to ensure the accuracy of your calcuation?

Footnotes

 [1] van Beest, B.W.H, Kramer, G.J., van Santen, R.A., PRL, 64, 16, 1990.
 [2] Harris, J. G., Yung, K. H., J. Phys. Chem. 1995, 99.